Stacking order, interaction and weak surface magnetism in layered graphene sheets 
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Recent transport experiments have demonstrated that the rhombohedral stacking trilayer 
graphene is an insulator with an intrinsic gap of 6meV and the Bernal stacking trilayer one is 
a metal. We propose a Hubbard model with a moderate U for layered graphene sheets, and show 
that the model well explains the experiments of the stacking dependent energy gap. The on-site 
Coulomb repulsion drives the metallic phase of the non-interacting system to a weak surface anti- 
ferromagnetic insulator for the rhombohedral stacking layers, but does not alter the metallic phase 
for the Bernal stacking layers. 
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In the past several years, the rapid development in 
preparing few layer graphene samples has promoted great 
theoretical [l-12| and experimental [l3-21| interests in 
such novel quasi-two-dimensional electron systems. The 
few layer graphene may be a platform for many new 
physics issues and is of potential application in electron- 
ics. One peculiar feature of the layered graphene system 
is the stacking order, which offers a new route to manip- 
ulate the electronic properties in graphene layers. 

The Bernal (or ABA) stacking and the rhombohedral 
(or ABC) stacking are two stable stacking orders ob- 
served in experiments. As shown in Fig. 1, in either 
ABA or ABC stacking order, the second graphene sheet 
is shifted by one bond length along the C-C bond direc- 
tion. The third graphene sheet is shifted back and aligned 
with the first sheet in the ABA stacking, while is shifted 
further by one more bond length along the same direc- 
tion in the rhombohedral stacking. So the ABA stacking 
order is ABABAB- • • , and the rhombohedral stacking is 
ABCABC- • • . The trilayer graphene system is the mini- 
mal structure relevant to the stacking orders. 

The electronic structures of the graphene layers 
strongly depend on their stacking orders [l|43|- In the 
ABA stacking TV-layer system, there are N/2 electron- 
like and N/2 hole-like parabolic sub-bands touching at 
e = for even N, and an additional sub-band with a lin- 
ear dispersion for odd N. The states in all the sub-bands 
are bulk states extended to all the layers. In the ABC 
stacking layers, the low energy electronic structure is de- 
scribed by two sub-bands with dispersion e ~ k N near 
the points K and K' in the 2D Brillouin zone. These 
low energy states are localized on the outermost layers, 
and are zero modes on the surfaces protected by the 
topology 0, [Hi- In two dimension, the dispersion of 
e(k) ~ k N gives a density of states D(e) ~ e~ 1 + 2 / iv ) 



which is divergent for N > 3 at e = 0. This indicates a 
strong instability toward symmetry broken states (3. [23l|. 
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FIG. 1: (Color online). Schematic diagrams of trilayer 
graphene sheets, (a): Bernal (ABA) stacking and (b): rhom- 
bohedral (ABC) stacking, (c) and (d) are their side views. 
Blue and pink colors represent carbon atoms on sub-lattices 
A and B, respectively. 

Trilayer graphene systems are of particular interest for 
it represents the simplest case for the stacking dependent 
graphene. Very recently, a stacking dependent intrinsic 
gap in trilayer graphene has been observed in the trans- 
port measurement [13, El- In the charge neutral case, 
namely undoped trilayer samples, the experiments indi- 
cate that the ABA stacking trilayer graphene is metallic, 
whereas the ABC stacking trilayer graphene is insulat- 
ing with an intrinsic gap about 6 meV. Since the non- 
interacting electronic structure of both stacking orders 
are gapless hence metallic, the experimental observation 
of the gap in ABC stacking trilayer is in sharp contrac- 
tion with the non-interaction picture and points to the 
importance of the interaction in these systems. 

In this Letter, we propose that the observed stacking- 
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dependent metallic or insulting states can be explained 
by a Hubbard model with a moderate on-site Coulomb 
repulsion U . We use a self-consistent mean field theory to 
show that the ground state of the ABC stacking trilayer 
is a weak anti-ferromagnet with opposite ferrimagnetic 
orderings on the top and bottom layers, due to the diver- 
gent density of states in the metallic phase. The magnetic 
ordering opens a gap, which is in good comparison with 
the experimental data. Our theory shows that the metal- 
lic phase of the ABA stacking trilayer is stable against 
a moderate Hubbard U due to the non-divergent den- 
sity of state. Our theory is extended to study stacking- 
dependent graphene systems for larger numbers of lay- 
ers. We have found that it is a general property for the 
ABC stacking graphene layers that the on-site Hubbard 
U opens a gap at the Fermi level and leads to a weak sur- 
face antifcrromagnetic state. Our results can be further 
tested in future experiments. 

We model N layer graphene systems by using a Hub- 
bard model H — Hq + Hu, where Hq = H lntrll + Hi ntcr is 
a tight binding Hamiltonian to describe the kinetic term 
of the system and Hjj describes the on-site Coulomb re- 
pulsion. The chemical potential is set to zero, and the 
average electron per site is one. The intralayer hopping 
term if intra is the tight-binding Hamiltonian of indepen- 
dent graphene sheets. For simplicity, we only include 
nearest neighbor hoppings [24| 



Y.WiMh<,(j) + h.c} (1) 

l(ij)cr 



where ai a (i) and bi a (j) are the annihilation operators 
of an electron on sublatticcs A and B, respectively. I 
denotes the layer index running from 1 to N, and (ij) 
nearest neighbor pairs, and a the spin. ifintcr describes 
the interlayer hopping given by 
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{al(i)b Va (i') + h.c}. (2) 



for the rhombohedral or Bernal stacking orders. Here, 
{IV) is summed over the two adjacent layers, and {W) is 
summed over two sites aligned in adjacent layers as shown 
in Fig. 1. The Hubbard term Hjj — U J2n n /t(*) n u(*) 
will be approximated by a mean field Hamiltonian, 
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where a = —a. {ni a {i)) is determined self-consistcntly. 
We consider a spin density wave state and introduce two 
mean fields on each layer I, one for sublattice A and one 
for sublattice B, {n^ B ). The mean fields for spin down 

are related to the spin-up ones, {nf^ B ) = 1 — {nf^ B ). 
Note that we have examined possible charge density wave 
states within the model and found no evidence for that. 

We first examine the trilayer graphenes (N = 3). The 
energy bands for the non-interacting models are shown 
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FIG. 2: (Color online). Low energy electronic bands in ABC 
(panel a) and ABA stacking (panel b) trilayer graphene at 
(7 = 0. Panel (c): mean field energy gaps as functions of 
U for ABC and ABA trilayer graphene. Panel (d): energy 
bands of ABC stacking with U = 6.2eV. The parameters are 
t = 3.16eV and tx = 0.39eV. 



in Fig. 2 (a) for ABC stacking and in Fig. (b) for ABA 
stacking with parameters t = 3.16eV and t± = 0.39eV. 
The dispersion in (a) is e ~ fc 3 at small k for both conduc- 
tion and valence bands, which gives rise to a divergent 
density of states D{t) — e -1 / 3 at e = 0, and the wave 
functions for k near the K or K' points are localized on 
the outer surfaces. The energy bands in (b) consist of a 
parabolic and a linear dispersions, both of which are not 
localized on the outer surfaces, and the density of states 
is a constant at e = 0. In the presence of the Hubbard U, 
the spin density wave ordering occurs at any U > for 
the ABC stacking case, and only at U > £/ c ABA « 6.3eV. 
The energy gaps associated with the spin density wave 
orderings are plotted in Fig. 2(c) as functions of U for 
both the ABC and ABA stacking orders. At U < C/ C ABA , 
the energy bands of the ABA stacking trilayer graphene 
are rigid against the Hubbard U. The sharp distinction 
between the ABC and ABA stacking trilayer graphene 
is attributed to their different density of states near the 
Fermi level. The divergent density of states of the surface 
zero modes for the ABC stacking graphene, protected by 
the momentum topology, actually induces its sensitivity 
to the interaction. 

More quantitatively, there are three distinguished re- 
gions in U for the gaps. At U < 5.5eV, the energy gap is 
zero for the ABA stacking and is tiny for the ABC stack- 
ing. At 5.5eV < U < 6.4eV, the gap size grows rapidly 
to be observable (several meV) for the ABC stacking, 
but remains zero or tiny for the ABA stacking. In this 
region, the ABC stacking trilayer is insulating with an 
observable gap while the ABA stacking trilayer remains 
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conduct. At U > 6.4eV, the gaps for both ABC and 
ABA stacking orders become observable, and become in- 
sulating. Actually the gaps for the two stacking orders 
become similar at U > 7eV as we can see from Fig. 2(c). 
Experimentally, the transport data shows ABC stacking 
trilayer graphene is an insulator with a gap of 6meV and 
the ABA stacking trilayer is metallic. In comparison with 
the experiments, the mean field calculations of the Hub- 
bard model suggest that the Hubbard U is within the 
interval of a moderate values 5.5eV < U < 6.4eV, i.e. 
1.74* < U < 2.03i 

In Fig. 2(d), we show the calculated quasi-particle 
dispersion for the ABC stacking trilayer graphene for 
a choice of U = 6.2eV. The corresponding gap is esti- 
mated to be e g ss 5.8meV. Our model and the calcu- 
lations well explain the recent experiments showing the 
stacking-dependent energy gap in trilayer graphene. The 
experimentally observed energy gap may be used to esti- 
mate the value of U. Our mean field theory suggests that 
U ks 6.2eV. More accurate numerical simulation may im- 
prove this estimate. 
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FIG. 3: (Color online). Schematic illustration of the spin or- 
derings (panel a) and calculated net spin polarization per site 
in each layer L (panel b) in the charge neutral ABC stacking 
trilayer graphene. In the calculations, U = 6.2eV, and the 
hopping parameters are the same as in Fig. 2. 

We now discuss the spin density state and the spin 
polarization of the ABC stacking layer. From the self- 
consistent mean field theory we obtain the site spin po- 
larization on sublattice A or B, defined as P a (l,i) = 



A,B 



) — {nf{ B ). The calculated spin polarizations arc 
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plotted in Fig. 3(b), and the spin structure in the tri- 
layer graphene is schematically illustrated in Fig. 3(a). 
The spin ordering is antiferromagnetic, where the neigh- 
boring spins (intra- or inter- layer) are anti-parallel to 
each other. However, there is a net spin polarization 
on the top or bottom layer, so each surface shows ferri- 
magnetic ordering. The spin polarization is mainly dis- 
tributed on the two outer surfaces. In each layer, the spin 
polarizations on sublattices A and B have opposite direc- 
tions. The net spin polarization is zero in the mid layer, 
and has opposite sign in the top or bottom layer. There 
is a symmetry of combined inversion and time reversal: 
P s (l = 1, i e A(B)) = P s (l = 3, i g B(A)). Note that the 
average spin polarization of the whole system is zero. For 
the parameters given in Fig. 3(b), the site spin polariza- 



tions in the top layer are about 6.9 x 10 -4 and —9.6 x 10 -4 
on sublattices A and B, respectively, and the net spin po- 
larization is —2.5 x 10~ 4 per site in average, which gives 
a surface magnetization 0.005/is/nm 2 . The weak sur- 
face magnetization on the ABC trilayer graphene is in 
analogy with the ferromagnetic edge states in graphene 
zigzag ribbon 25j, in which the density of states of the 
flat band edge states is divergent, inducing the edge spin 
polarization in the presence of a weak interaction. The 
interaction induced gap in graphene zigzag ribbon has 
been confirmed in a recent STM experiment 
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FIG. 4: (Color online). ABC stacking JV-layer graphene. 
Energy band (panel a) and spin polarization (panel b) ob- 
tained in mean field theory for N = 24. The gaps e B i and e 9 2 
(panel c) and spin polarization (panel d) as functions of N. 
U — Q.2eV , and the hopping parameters are the same as in 
Fig. 2. 

We now discuss N > 3 graphene layers. For the ABC 
stacking graphene layers and in the charge neutral case, 
there is always an interaction (U) induced gap at the 
Fermi level with spontaneous surface spin density wave 
ordering. In Fig. 4(a) and (b), we show the results of 
N = 24 layers as an example. Since the minimum band 
gap (e g i) is no longer at the K or K' points as we can see 
from Fig. 4(a), we introduce the second gap e g 2 for the 
energy gap at the K or K' points. As shown in Fig. 4(b), 
the spin polarization is localized near the surfaces. In Fig. 
4 (c), we present the energy gaps as functions of N. As 
N increases, e g i first increases to approach its maximum 
about 20meV at N — 9, then decreases to a value of 
15meV at N = 24. On the other hand, e g 2 increases with 
the layer thickness N and reaches a saturated value about 
33mcV. The iV-dependent spin polarization is shown in 
Fig. 4 (d), which increases with the layer thickness, and 
approaches to a saturated value, which is at least 5 to 
6 times of the surface magnetization in the trilayer case. 
Note that the first principle calculations involving a local 
spin density approximation has been applied to study 8 
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layers ABC stacking graphene [12j, and reported a spin 
density wave ground state. Their result is consistent with 
the results of the Hubbard model proposed here, while 
our results are more general and distinguishes different 
stacking orders. 

We have also applied the mean field theory to study 
iV-layer graphene of the ABA stacking. The results are 
similar to the trilayer case and there is a threshold U ABA , 
which is weakly TV dependent and remains to be finite at 
large N, the metallic phase is stable against the spin 
density wave ordering at U < U ABA . 

We argue that the mean field theory should give a 
qualitatively or semi-quantitatively correct physics on 
the stacking dependent instability, or the insulating or 
metallic states in layered graphene, while more accurate 
calculations may refine the estimate of the value of U. We 
remark that the proposed Hubbard model with a moder- 
ate U should capture the most important physics for the 
stacking dependent ground states in layered graphene. 
The intersite Coulomb repulsion has tendency to drive 
the metallic phase to a charge density wave state, which is 
not compatible with the on-site U studied in the present 
work. Since the intersite repulsion is relatively weaker 
than the on-site Coulomb repulsion U, we may argue that 
that term may not be relevant. More exotic states such 
as quantum spin Hall state and anomalous Hall state 
have been proposed in models with spin-orbit coupling 
or intersite interaction on honeycomb lattice [2^, HH. 
The possible realization of these exotic phases in layered 
graphene will be highly interesting. In view of the very 
weak spin-orbit coupling in graphene [29], more detailed 
study will be needed to explore the possibility. 

In summary, we have proposed a Hubbard model with 
a moderate U to describe N layer graphene, and applied 
a mean field theory to study the ground state and the 
excited energy gap of the charge neutral systems. The 
ground state of the non-interacting graphene layers are 
metallic, whose density of states at the Fermi level is 
divergent as D(e) ~ e - 1 + 2 / N for rhombohedral (ABC) 
stacking, and is a constant for Bernal (ABA) stacking. 
The metallic ground state of the ABA stacking layer is 
stable against the on-site Coulomb repulsion for moder- 
ate value of U below a threshold U ABA about 6-7 eV. 
The metallic state of the ABC stacking layer is found to 
be unstable against any repulsion U due to the divergent 
density of states at zero energy. Its ground state is surface 
antiferromagnetic state with opposite ferrimagnetism on 
top or bottom surfaces, which opens a gap. The energy 
gap is estimated to be 5.8meV for U = 6.2eV for N = 3. 
Our model and calculations well explain the recent trans- 
port experiments, showing that the ABC stacking trilayer 
graphene is an insulator with a gap about 6 meV, and 
ABA stacking trilayer graphene remains to be metallic. 
The spin polarization in the spin ordered state is found 
to be weak, but should be measurable. 
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